% no Gaussian modes in hostograms!

load('colors.mat');

%% para
% measDir = ['..' filesep '210518_Faruq_1Wildbat'];

% flagSaveFigs = 0;
% flagGate = 1;

% rSqrCutoff = 0.6;
% areaGateRange = [10 30];
% velocityGateRange = [2 5];
% tau1GateRange = [0 12];
% tau2GateRange = [0 12];
% deltaDefoRange = [-.1 .5];

%% find samples and choose
dirRes = dir([measDir filesep 'FTFitData' filesep]);
dirRes = dirRes([dirRes.isdir]);
samples = cell(1,length(dirRes)-2);
flowRates = cell(1,length(dirRes)-2);
fprintf('Index\tName\t\t\tFlowRate\n');
for k = 3:length(dirRes)
    idxUnder = find(dirRes(k).name == '_');
    samples{k-2} = dirRes(k).name(1:idxUnder(end-1)-1);
    flowRates{k-2} = dirRes(k).name(idxUnder(end-1)+1:idxUnder(end)-1);
    fprintf('%d\t%s\t\t%s\n',k-2,samples{k-2},flowRates{k-2});
end
% idx = input('Choose samples: ');
% if not(isempty(idx))
%     samples = samples(idx);
%     flowRates = flowRates(idx);
% end


%%
for k = 1:6
    figure(k);clf;
end

fprintf([repmat('-',1,length(samples)) '\n']);
for idx = 1:length(samples)
dirRes = dir([measDir filesep 'FTFitData' filesep samples{idx} '_' flowRates{idx} '_*' ]);

rootDir = [dirRes.folder filesep dirRes.name];
load([rootDir filesep 'data.mat']);
idxs = find(filedir == filesep);
sampleName = filedir(idxs(end-2)+1:idxs(end-1)-1);

%% area & velocity vectors
areaVec = zeros(size(cellNums));
veloVec = zeros(size(cellNums));
ptr = 1;
for k = 1:length(sampleData)
    tmp = [sampleData(k).measData.cellData.data];
    if ~isempty(tmp)
        areaVec(ptr:ptr+length(tmp)-1) = [tmp.medianArea];
        veloVec(ptr:ptr+length(tmp)-1) = [tmp.medianVelocity];
        ptr = ptr + length(tmp);
    end
end

%% deltaDefo
deltaDefoTau1 = defoTau1Inlet - defoTau1Outlet;
deltaDefoTau2 = defoTau2Outlet - defoTau2Inlet;

%% gates
maskTau1Consistency = (taus1Circularity >= 0) & (rSqrsCircularityTau1 >= rSqrCutoff);
maskTau2Consistency = (taus2Circularity >= 0) & (rSqrsCircularityTau2 >= rSqrCutoff);

if flagGate
    maskArea = (areaGateRange(1) <= areaVec) & (areaVec <= areaGateRange(2));
    maskVelo = (velocityGateRange(1) <= veloVec) & (veloVec <= velocityGateRange(2));
    maskCell = maskArea & maskVelo;
    
    maskTau1Range = (tau1GateRange(1) <= taus1Circularity) & (taus1Circularity <= tau1GateRange(2));
    maskTau2Range = (tau2GateRange(1) <= taus2Circularity) & (taus2Circularity <= tau2GateRange(2));
    maskDeltaDefoTau1Range = (deltaDefoRange(1) <= deltaDefoTau1) & (deltaDefoTau1 <= deltaDefoRange(2));
    maskDeltaDefoTau2Range = (deltaDefoRange(1) <= deltaDefoTau2) & (deltaDefoTau2 <= deltaDefoRange(2));
    
    strGate = '_gate';
else
    maskCell = true(size(areaVec));
    
    maskTau1Range = true(size(taus1Circularity));
    maskTau2Range = true(size(taus1Circularity));
    maskDeltaDefoTau1Range = true(size(deltaDefoTau1));
    maskDeltaDefoTau2Range = true(size(deltaDefoTau2));
    
    strGate = '_noGate';
end

recFig = 0;

%% area
recFig = recFig + 1; figure(recFig);clf;
histogram(areaVec(maskCell));
xlabel('Area (um^2)');
title([sampleName ' @ ' flowRate ': ' sprintf('%d cells',sum(maskCell))],'Interpreter','none');
if flagSaveFigs
    saveas(recFig,[rootDir filesep 'Area' strGate '.fig']);
end
print(recFig, '-dpng','-painters',[rootDir filesep 'Area' strGate]);

%% velocity
recFig = recFig + 1; figure(recFig);clf;
histogram(veloVec(maskCell));
xlabel('Velocity (cm/s)');
title([sampleName ' @ ' flowRate ': ' sprintf('%d cells',sum(maskCell))],'Interpreter','none');
if flagSaveFigs
    saveas(recFig,[rootDir filesep 'Velocity' strGate '.fig']);
end
print(recFig, '-dpng','-painters',[rootDir filesep 'Velocity' strGate]);

%% tau1
recFig = recFig + 1; figure(recFig);clf;
subplot(4,1,1);
edges = 0:0.1:1;
histogram(rSqrsCircularityTau1(maskCell & maskTau1Consistency), 'BinEdges',edges);
hold on;
histogram(rSqrsCircularityTau1(maskCell & maskTau1Consistency & maskTau1Range), 'BinEdges',edges, 'facealpha',.5);
hold off;
yLim = get(gca,'YLim');
hold on; line(rSqrCutoff*[1,1],yLim,'Color',color1); hold off;
set(gca,'YLim',yLim);
ylabel('# of R^2');
title([sampleName ' @ ' flowRate ': ' sprintf('%d cells',sum(maskCell & maskTau1Consistency & maskTau1Range))],'Interpreter','none');

subplot(4,1,2:4);
histogram(log10(taus1Circularity(maskCell & maskTau1Consistency & maskTau1Range)));
xlabel('\tau_1 (ms)');
ylabel('#');
logXAxisLabel();

if flagSaveFigs
    saveas(recFig,[rootDir filesep 'Tau1' strGate '.fig']);
end
print(recFig, '-dpng','-painters',[rootDir filesep 'Tau1' strGate]);

%% deltaDefoTau1
recFig = recFig + 1; figure(recFig);clf;
histogram(deltaDefoTau1(maskCell & maskTau1Consistency & maskTau1Range & maskDeltaDefoTau1Range));
xlabel('\Deltadefo(\tau_1)');
ylabel('#');
title([sampleName ' @ ' flowRate ': ' sprintf('%d cells',sum(maskCell & maskTau1Consistency & maskTau1Range & maskDeltaDefoTau1Range))],'Interpreter','none');

if flagSaveFigs
    saveas(recFig,[rootDir filesep 'DeltaDefoTau1' strGate '.fig']);
end
print(recFig, '-dpng','-painters',[rootDir filesep 'DeltaDefoTau1' strGate]);

%% tau2
recFig = recFig + 1; figure(recFig);clf;
subplot(4,1,1);
edges = 0:0.1:1;
histogram(rSqrsCircularityTau2(maskCell & maskTau2Consistency), 'BinEdges',edges);
hold on;
histogram(rSqrsCircularityTau2(maskCell & maskTau2Consistency & maskTau2Range), 'BinEdges',edges, 'facealpha',.5);
hold off;
yLim = get(gca,'YLim');
hold on; line(rSqrCutoff*[1,1],yLim,'Color',color1); hold off;
set(gca,'YLim',yLim);
ylabel('# of R^2');
title([sampleName ' @ ' flowRate ': ' sprintf('%d cells',sum(maskCell & maskTau2Consistency & maskTau2Range))],'Interpreter','none');

subplot(4,1,2:4);
histogram(log10(taus2Circularity(maskCell & maskTau2Consistency & maskTau2Range)));
xlabel('\tau_2 (ms)');
ylabel('#');
logXAxisLabel();

if flagSaveFigs
    saveas(recFig,[rootDir filesep 'Tau2' strGate '.fig']);
end
print(recFig, '-dpng','-painters',[rootDir filesep 'Tau2' strGate]);

%% deltaDefoTau2
recFig = recFig + 1; figure(recFig);clf;
histogram(deltaDefoTau2(maskCell & maskTau2Consistency & maskTau2Range & maskDeltaDefoTau2Range));
xlabel('\Deltadefo(\tau_2)');
ylabel('#');
title([sampleName ' @ ' flowRate ': ' sprintf('%d cells',sum(maskCell & maskTau2Consistency & maskTau2Range & maskDeltaDefoTau2Range))],'Interpreter','none');

if flagSaveFigs
    saveas(recFig,[rootDir filesep 'DeltaDefoTau2' strGate '.fig']);
end
print(recFig, '-dpng','-painters',[rootDir filesep 'DeltaDefoTau2' strGate]);


fprintf('-');
end
fprintf('\n');

return
%% functions
function logXAxisLabel()
range = get(gca,'XLim');
xTickVals = [1 2 5]' * 10.^(floor(range(1)):ceil(range(2)));
xTickVals = xTickVals(:);
xTickVals = xTickVals( (10.^range(1) <= xTickVals) & (xTickVals <= 10.^range(2)) );
set(gca,'XTick',log10(xTickVals));
set(gca,'XTickLabel',cellstr(num2str(xTickVals, '%-g')));

end